
#include <vector>
#include <cmath>
#include <limits>
#include <functional>
#include <algorithm>

#include "member_BSplineCurve.h"

using namespace BWE;

member_BSplineCurve::member_BSplineCurve(BSplineCurve* polyStripe)
{
	boss = polyStripe;
	thickness = 0;
	closeLoop = false;
	dirty = false;
	length = 0;
}
member_BSplineCurve::~member_BSplineCurve()
{

}

inline void thomasSolve(BArray<BReal>& a, BArray<BReal>& b, BArray<BReal>& c, BArray<BReal>& d, BArray<BReal>& x)
{
    int n = a.size();
    BArray<BReal> cp(c);
    BArray<BReal> dp(d);

    for (int i = 1; i < n; ++i)
    {
        BReal m = a[i] / b[i - 1];
        b[i] -= m * cp[i - 1];
        dp[i] -= m * dp[i - 1];
    }

    x[n - 1] = dp[n - 1] / b[n - 1];
    for (int i = n - 2; i >= 0; --i)
    {
        x[i] = (dp[i] - cp[i] * x[i + 1]) / b[i];
    }
}

BArray<BReal> member_BSplineCurve::computeCoefficients(const BArray<BReal>& values)
{
    const int n = values.size();
    BArray<BReal> result(n);
    if (n <= 2)
        return result;

    BArray<BReal> rhs(n);
    for (int i = 1; i < n - 1; ++i)
    {
        rhs[i - 1] = 6 * (values[i + 1] - 2 * values[i] + values[i - 1]);
    }

    BArray<BReal> a(BReal(4.0), n);
    BArray<BReal> b(BReal(1.0), n);
    BArray<BReal> c(BReal(1.0), n);

    for (int i = 1; i < n - 3; ++i)
    {
        BReal factor = c[i - 1] / a[i - 1];
        a[i] -= factor * b[i - 1];
        rhs[i] -= factor * rhs[i - 1];
    }

    result[n - 2] = rhs[n - 3] / a[n - 3];
    for (int i = n - 4; i >= 0; --i)
    {
        result[i + 1] = (rhs[i] - b[i] * result[i + 2]) / a[i];
    }
    return result;
}
